Max criteria above 250 tons

fac <- mutate(fac, 
              max_value = max_crit,
              max_crit  = ifelse(max_crit > 249.9, 265, max_crit))

fac <- mutate(fac, max_crit = ifelse(max_crit < 5, 5, max_crit))

fac_crit <- group_by(fac, SOURCE_ID) %>% 
            summarize(max_crit = max_crit[1]) %>%
            arrange(max_crit) %>%
            mutate(id = 1:n())

#print(nrow(filter(fac_crit, max_crit > 249.9)))

fac_crit$colors <- NA

for(i in 1:nrow(fac_crit)) fac_crit$colors[i] <- (fac_crit$max_crit[i] > 249.9) * runif(1) > 0.6


ggplot(fac_crit, aes(id, max_crit)) + 
  geom_bar(stat="identity", aes(fill = max_crit > 249.9)) + 
  scale_fill_manual(values = ton_colors, 
                    name   = "",
                    labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(fill = guide_legend(ncol = 2, keywidth=3)) +
  theme(legend.position="top")

*120 facilities out of 2,000

ggplot(fac_crit, aes(id, max_crit)) + 
       geom_bar(stat="identity", aes(fill = colors)) + 
       scale_fill_manual(values= ton_colors, 
                         name="",
                         labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(fill = guide_legend(ncol = 2, keywidth=3)) +
  theme(legend.position="top")

Max criteria vs. total HAPs

fac  <- group_by(fac, SOURCE_ID, Facility) %>% 
        summarize(max_crit     = max_crit[1],
                  max_value    = max_value[1],
                  total_HAPS   = sum(Emissions, na.rm=T),
                  Formaldehyde = sum(Emissions[Pollutant == "FORMALDEHYDE"], na.rm=T),
                  Acrolein     = sum(Emissions[Pollutant == "ACROLEIN"], na.rm=T),
                  Lead         = sum(Emissions[Pollutant == "LEAD"], na.rm=T),
                  Nickel       = sum(Emissions[Pollutant == "NICKEL"], na.rm=T),
                  Arsenic      = sum(Emissions[Pollutant == "ARSENIC"], na.rm=T),
                  Metals       = sum(Emissions[Pollutant %in% c("LEAD", "NICKEL", "ARSENIC")], na.rm=T))

fac <- filter(fac, max_value > 1, total_HAPS > 0.05) %>%
       ungroup() %>%
       arrange(max_value) 

fac$id <- 1:nrow(fac)

fac$colors <- as.factor(fac$max_value > 249.9)


# HAPs arranged in order of increasing Criteria
ggplot(fac, aes(id, log(total_HAPS))) + geom_point(size=2.8 * size_X, alpha=0.65) + labs(title = "Facility HAPs emissions in order of increasing criteria emissions")

# Total HAPs
p <- ggplot(fac) +
     labs(title = "Total HAPs") +
     geom_point(data = fac, aes(id, log10(total_HAPS), color = colors), size = 2.5 * size_X, alpha = 0.7) +
     scale_color_manual(values= ton_colors, 
                         name="",
                         labels = c(" Below 250 tons  ", " Above 250  ")) +
    guides(color= guide_legend(ncol = 2, keywidth=3)) +
    theme(legend.position="top")


ggplot(fac) + 
  labs(title = "HAPs emissions vs maximum Criteria emissions") +
  geom_point(data = fac, aes(log10(max_value), log10(total_HAPS), color = colors), size = 3.5 * size_X, alpha = 0.7) +
  geom_smooth(aes(log10(max_value), log10(total_HAPS)), span = 155) + #method="lm")
  scale_color_manual(values= ton_colors, 
                         name="",
                         labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Top 10% HAPS emitters

ggplot(filter(fac, total_HAPS > quantile(total_HAPS, 0.9)[[1]])) + 
  labs(title = "") +
  geom_point(data = filter(fac, total_HAPS > quantile(total_HAPS, 0.9)[[1]]), aes(log10(max_value), log10(total_HAPS), color = colors), size = 4 * size_X, alpha = 0.7) +
  scale_color_manual(values= ton_colors, 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Formaldehyde

ggplot(filter(fac, Formaldehyde > 0.00001)) + 
  labs(title = "Formaldehyde emissions vs maximum Criteria emissions") +
  geom_point(data = filter(fac, Formaldehyde > 0.00001), aes(log10(max_value), log10(Formaldehyde), color = colors), size = 3.5 * size_X, alpha = 0.7) +
  geom_smooth(aes(log10(max_value), log10(Formaldehyde)), span = 15) + 
  scale_color_manual(values= ton_colors, 
                         name="",
                         labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Top 10% formaldehyde emitters

ggplot(filter(fac, Formaldehyde > quantile(Formaldehyde, 0.9)[[1]])) + 
  labs(title = "") +
  geom_point(data = filter(fac, Formaldehyde > quantile(Formaldehyde, 0.9)[[1]]), aes(log10(max_value), log10(Formaldehyde), color = colors), size = 4 * size_X, alpha = 0.7) +
  scale_color_manual(values= ton_colors, 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Acrolein

ggplot(filter(fac, Acrolein > 1E-5)) + 
  labs(title = "Acrolein emissions vs maximum Criteria emissions") +
  geom_point(data = filter(fac, Acrolein > 1E-5), 
             aes(log10(max_value), log10(Acrolein), color = colors), size = 3.5 * size_X, alpha = 0.7) +
  geom_smooth(aes(log10(max_value), log10(Acrolein)), span = 15) + 
  scale_color_manual(values= ton_colors, 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Top 10% acrolein emitters

ggplot(filter(fac, Acrolein > quantile(Acrolein, 0.9)[[1]])) + 
  labs(title = "") +
  geom_point(data = filter(fac, Acrolein > quantile(Acrolein, 0.9)[[1]]), aes(log10(max_value), log10(Acrolein), color = colors), size = 4 *size_X, alpha = 0.7) +
  scale_color_manual(values= ton_colors, 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Lead

ggplot(filter(fac, Lead > 1E-5)) + 
  labs(title = "Lead emissions vs maximum Criteria emissions") +
  geom_point(data = filter(fac, Lead> 1E-5), 
             aes(log10(max_value), log10(Lead), color = colors), size = 3.5 * size_X, alpha = 0.7) +
  geom_smooth(aes(log10(max_value), log10(Lead)), span = 15) + 
  scale_color_manual(values= ton_colors, 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Top 10% lead emitters

ggplot(filter(fac, Lead > quantile(Lead, 0.9)[[1]])) + 
  labs(title = "") +
  geom_point(data = filter(fac, Lead > quantile(Lead, 0.9)[[1]]), aes(log10(max_value), log10(Lead), color = colors), size = 4 * size_X, alpha = 0.7) + 
  scale_color_manual(values= ton_colors, 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Nickel

ggplot(filter(fac, Nickel > 1E-5)) + 
  labs(title = "Nickel emissions vs maximum Criteria emissions") +
  geom_point(data = filter(fac, Nickel> 1E-5), 
             aes(log10(max_value), log10(Nickel), color = colors), size = 3.5 * size_X, alpha = 0.7) +
  geom_smooth(aes(log10(max_value), log10(Nickel)), span = 15) + 
  scale_color_manual(values= ton_colors, 
                         name="",
                         labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Top 10% nickel emitters

ggplot(filter(fac, Nickel > quantile(Nickel, 0.9)[[1]])) + 
  labs(title = "") +
  geom_point(data = filter(fac, Nickel > quantile(Nickel, 0.9)[[1]]), aes(log10(max_value), log10(Nickel), color = colors), size = 4.5 * size_X, alpha = 0.7) +
  scale_color_manual(values= c(viridis(10)[8], magma(5)[4]), 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top") 

Arsenic

ggplot(filter(fac, Arsenic > 1E-5)) + 
  labs(title = "Arsenic emissions vs maximum Criteria emissions") +
  geom_point(data = filter(fac, Arsenic > 1E-5), 
             aes(log10(max_value), log10(Arsenic), color = colors), size = 3.5 * size_X, alpha = 0.7) +
  geom_smooth(aes(log10(max_value), log10(Arsenic)), span = 15) + 
  scale_color_manual(values= ton_colors, 
                     name="",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Top 10% arsenic emitters

ggplot(filter(fac, Arsenic > quantile(Arsenic, 0.9)[[1]])) + 
  labs(title = "") +
  geom_point(data = filter(fac, Arsenic > quantile(Arsenic, 0.9)[[1]]), aes(log10(max_value), log10(Arsenic), color = colors), size = 4 * size_X, alpha = 0.7) +
  scale_color_manual(values = ton_colors, 
                     name   = "",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Metals

ggplot(filter(fac, Metals > 1E-4)) + 
    labs(title = "Metals emissions vs maximum Criteria emissions") +
    geom_point(data = filter(fac, Metals > 1E-4), aes(log10(max_value), log10(Metals), color = colors), size = 3.5 * size_X, alpha = 0.7) +
    geom_smooth(aes(log10(max_value), log10(Metals)), span = 15) + 
    scale_color_manual(values= ton_colors, 
                       name="",
                       labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Top 10% Metal emitters

ggplot(filter(fac, Metals > quantile(Metals, 0.9)[[1]])) + 
  labs(title = "") +
  geom_point(data = filter(fac, Metals > quantile(Metals, 0.9)[[1]]), aes(log10(max_value), log10(Metals), color = colors), size = 4 * size_X, alpha = 0.7) +
  scale_color_manual(values = ton_colors, 
                     name   = "",
                     labels = c(" Below 250 tons  ", " Above 250  ")) +
  guides(color = guide_legend(ncol = 2, keywidth=3, keyheight = 3, override.aes = list(size=8))) +
  theme(legend.position="top")

Table


library(knitr)

options(digits = 2)

fac[ ,c(4:11)] <- round(fac[ ,c(4:11)], 3)

DT::datatable(fac[ ,-c(3,12:13)], options=list(searching=F, paging=T, scrollX=T), rownames = FALSE)
#kable(fac[ ,-c(3,12:13)])